{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "provenance": []
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "source": [
        "![Automekin.png]()\n",
        "# Automated reaction Mechanisms and Kinetics\n",
        "\n"
      ],
      "metadata": {
        "id": "vXtSvTouM1dG"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Installing third-party packages"
      ],
      "metadata": {
        "id": "9LmhQ0GeMi9z"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "%%capture\n",
        "!pip install ase pubchempy\n",
        "!apt update\n",
        "!apt install bc gawk gcc gfortran parallel sqlite3 cm-super dvipng texlive-latex-extra texlive-latex-recommended\n",
        "%cd /opt/\n",
        "!git clone https://github.com/dgarayr/amk_tools.git\n",
        "%cd /opt/amk_tools\n",
        "!pip install -e ."
      ],
      "metadata": {
        "id": "cQDTrDFQGe3p"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Installing AutoMeKin"
      ],
      "metadata": {
        "id": "qTg3nK8uNBMr"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "P6BRAe0KGS70"
      },
      "outputs": [],
      "source": [
        "%%capture\n",
        "%cd /content/\n",
        "!git clone https://github.com/emartineznunez/AutoMeKin.git\n",
        "%cd /content/AutoMeKin\n",
        "! autoreconf -i\n",
        "!./configure --prefix=/opt/AutoMeKin\n",
        "!make && make install\n",
        "#We make use of the total number of processors\n",
        "!sed -i 's@ignore=1@ignore=0@g' /opt/AutoMeKin/bin/utils.sh\n",
        "%env PATH=\".:/opt/amk_tools/scripts:/opt/AutoMeKin/bin:/opt/AutoMeKin/bin/HLscripts:/opt/AutoMeKin/bin/MOPAC_DEV:/opt/bin:/usr/local/nvidia/bin:/usr/local/cuda/bin:/usr/local/sbin:/usr/local/bin:/usr/sbin:/usr/bin:/sbin:/bin:/tools/node/bin:/tools/google-cloud-sdk/bin\"\n",
        "%env LIBRARY_PATH=\"/opt/AutoMeKin/lib:/opt/AutoMeKin/bin/MOPAC_DEV:/usr/local/cuda/lib64/stubs\"\n",
        "%env AMK=/opt/AutoMeKin\n",
        "%env inter=0\n",
        "%cd /content\n",
        "!rm -rf /content/AutoMeKin"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Making the input files"
      ],
      "metadata": {
        "id": "13ZPJzKRQACQ"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "In this example the default molecule is **formic acid (HCOOH)** and it takes a few minutes (<8). **If you want to change the molecule**, you can just pick one (by **cid**, **name**, **smiles**, **sdf**, **inchi**, **inchikey** or **formula**) from: https://pubchem.ncbi.nlm.nih.gov/"
      ],
      "metadata": {
        "id": "iyHjiBHnUBHz"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "#@title Choose the default, formic acid (FA) molecule, or a new system.\n",
        "%%capture\n",
        "\n",
        "%cd /content\n",
        "!rm -rf test && mkdir test\n",
        "%cd /content/test\n",
        "\n",
        "import pubchempy as pcp\n",
        "\n",
        "#Building the molecule\n",
        "\n",
        "#The default is formic acid, but you might want another one https://pubchem.ncbi.nlm.nih.gov/\n",
        "sel = input('This example is for \\'formic acid\\'.Do you want to choose your own molecule (y/n)?: ')\n",
        "if sel == 'n':\n",
        "  molecule = 'formic acid'\n",
        "  nm = 'name'\n",
        "else:\n",
        "  nm = input('Identifier type: cid, name, smiles, sdf, inchi, inchikey, or formula: ')\n",
        "  molecule = input('Type your molecule here: ')\n",
        "\n",
        "query = pcp.get_compounds(molecule,nm,record_type='3d')\n",
        "m = query[0]\n",
        "f = open('mol.xyz','w')\n",
        "f.write(str(len(m.atoms)) + '\\n\\n')\n",
        "for i,a in enumerate(m.atoms): f.write(m.elements[i]+' '+str(a.x)+' '+str(a.y)+' '+str(a.z)+'\\n')\n",
        "f.close()\n",
        "\n",
        "#Fetching input file template\n",
        "!curl -L https://github.com/emartineznunez/AutoMeKin/raw/main/examples/FA.dat -o mol.dat\n",
        "!sed -i 's@FA@mol@;s@imagmin 200@imagmin 50@' mol.dat"
      ],
      "metadata": {
        "id": "2tOyyBzmQFYZ"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Running the test\n",
        "Running 3 iterations of the workflow with 5 tasks per iteration. Details of the program execution, input and output files, etc., can be looked up here: https://emartineznunez.github.io/AutoMeKin/docs/tutorial.html"
      ],
      "metadata": {
        "id": "FoE0ellxNoC-"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "!llcalcs.sh mol.dat 5 3 2"
      ],
      "metadata": {
        "id": "cZ2ysrp0KrSL"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "*Transition-state structures found in these 3 iterations*"
      ],
      "metadata": {
        "id": "1U8s9Tiy1X51"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "!tsll_view.sh"
      ],
      "metadata": {
        "id": "Cq2FHWHVZfrI"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "*Convergence table*"
      ],
      "metadata": {
        "id": "mpX1rSQN1qUz"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "!cat FINAL_LL_mol/convergence.txt"
      ],
      "metadata": {
        "id": "eCdAHDeccz3W"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "*Coarse-grained reaction network*"
      ],
      "metadata": {
        "id": "TXLRN_kW1jCY"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "!cat FINAL_LL_mol/RXNet.cg"
      ],
      "metadata": {
        "id": "i3iTpYbL1w7e"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "*Transition states and minima*"
      ],
      "metadata": {
        "id": "DNb_HAkC246Q"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "!awk 'BEGINFILE {print \"\\n\"}{print}' FINAL_LL_mol/TSinfo FINAL_LL_mol/MINinfo"
      ],
      "metadata": {
        "id": "VGc3Uy352gpa"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "*Kinetics: branching ratios and populations*"
      ],
      "metadata": {
        "id": "L6ou49Le2w3b"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "#Printing branching ratios and reading the csv file to plot the population vs time\n",
        "%cd /content/test\n",
        "#!awk '{print $0};NF==0{exit}' FINAL_LL_mol/kineticsE150\n",
        "#!awk '/Time/,0 {for(i=1;i<=NF-1;i++) printf \"%s, \",$i;print $NF}' FINAL_LL_mol/kineticsE150 > FINAL_LL_mol/kinetics.csv\n",
        "\n",
        "# importing libraries and creating the plot\n",
        "import pandas as pd\n",
        "from matplotlib import pyplot\n",
        "\n",
        "pyplot.rcParams['text.usetex'] = True\n",
        "pyplot.xticks(fontsize=14)\n",
        "pyplot.yticks(fontsize=14)\n",
        "\n",
        "data = pd.read_csv('FINAL_LL_mol/kinetics.csv')\n",
        "for count,col in enumerate(data.columns):\n",
        "  if count == 0:\n",
        "    x=data[col]\n",
        "    xcol=col\n",
        "  if count >=1: pyplot.plot(x,data[col],label=col,linewidth=1.0)\n",
        "\n",
        "pyplot.ylabel('Population',fontsize=20)\n",
        "pyplot.xlabel(xcol,fontsize=20)\n",
        "\n",
        "pyplot.legend()\n",
        "pyplot.xlim(0,max(x))\n",
        "pyplot.ylim(0,1000)\n",
        "pyplot.tight_layout()\n",
        "pyplot.savefig('/content/test/kinetics.png')\n",
        "pyplot.show()"
      ],
      "metadata": {
        "id": "XQJqcs5T2fZl"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "*Generate network.html for visualization*"
      ],
      "metadata": {
        "id": "fX1nuW9k13DS"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "from os import path\n",
        "\n",
        "%cd /content/test\n",
        "if path.isfile('tsdirLL_mol/KMC/starting_minimum'): tag = !awk 'NR=1{print $1}' tsdirLL_mol/KMC/starting_minimum\n",
        "else: tag = [1]\n",
        "!amk_gen_view.py FINAL_LL_mol RXNet.cg --b --paths MIN{tag[0]}\n",
        "!mkdir HTML\n",
        "!mv network.html HTML\n",
        "%cd /content/test/HTML"
      ],
      "metadata": {
        "id": "pTH02rKqmGSm"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# @title ##### _Open the file network.html_\n",
        "import IPython\n",
        "from IPython.core.magic import register_line_magic\n",
        "import subprocess\n",
        "\n",
        "@register_line_magic\n",
        "def run_local_server(line):\n",
        "    handle = IPython.display.display(\n",
        "            IPython.display.Pretty(\"Launching my server...\"),\n",
        "            display_id=True,\n",
        "    )\n",
        "    subprocess.Popen(['python', '-m', 'http.server'])\n",
        "    shell = \"\"\"\n",
        "        (async () => {\n",
        "            const url = new URL(await google.colab.kernel.proxyPort(8000, {'cache': true}));\n",
        "            const iframe = document.createElement('iframe');\n",
        "            iframe.src = url;\n",
        "            iframe.setAttribute('width', '100%');\n",
        "            iframe.setAttribute('height', '800');\n",
        "            iframe.setAttribute('frameborder', 0);\n",
        "            document.body.appendChild(iframe);\n",
        "        })();\n",
        "    \"\"\"\n",
        "    script = IPython.display.Javascript(shell)\n",
        "    handle.update(script)\n",
        "!fuser -k 8000/tcp\n",
        "%run_local_server"
      ],
      "metadata": {
        "cellView": "form",
        "id": "Hk9n87rhkXZJ"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "*Statistics of the reaction network*"
      ],
      "metadata": {
        "id": "Y6czZObNQbCP"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "%cd /content/test\n",
        "!amk_rxn_stats.py FINAL_LL_mol\n",
        "%cat rxn_stats.txt"
      ],
      "metadata": {
        "id": "RL01GCjJP-XG"
      },
      "execution_count": null,
      "outputs": []
    }
  ]
}